#include "reconstruction.hpp"

namespace hms::swe
{

template<typename MeshType>
template<BlockLocation loc>
void SolverThread<MeshType>::computeFluxes(
	const Block<loc>& block,
	const Array3XsMap& q,
	      Array3XsMap& fluxes,
	const Array1XsMap& h,
	const Ref<const ArrayXXs>& d,
	const Ref<const ArrayXXs>& z,
	const Ref<const ArrayXXs>& v,
	const MeshBase<MeshType>& mesh
){
	using RT = ReconstructionType;
	using ER = EdgeRange;

	// std::stringstream s;
	// s << "Thread #" << omp_get_thread_num() << ", edges ";

	#define HLLC(REC_TYPE, EDGE_RANGE, EDGE_INDEX, LR_BLOCK) \
	fluxScheme<REC_TYPE, EDGE_RANGE>( \
		block, q, fluxes, h, d, z, v, mesh, EDGE_INDEX, LR_BLOCK)
		// block, std::forward<Ts>(args) ..., EDGE_INDEX, LR_BLOCK)
	
		// ; ++(this->nEdges)
		// ; s << EDGE_INDEX << " "
		// << ( EDGE_RANGE == ER::horizontal ? 'h' : 'v' )
		// << "(" << (REC_TYPE == RT::secondOrder ? 2 : 1) << ". o.), lr: " << (LR_BLOCK).transpose() << "\n"

	/* the inner edges loop goes from x0-nx and y0-ny. If the block location is
	 * not "inner", then these must be modified and additional loops run. */
	constexpr Index
		x0 { atRight<loc>() ? 1 : 0 },
		y0 { atTop  <loc>() ? 1 : 0 };
	Index
		nx { block.nx() - ( atLeft  <loc>() ? 1 : 0 ) },
		ny { block.ny() - ( atBottom<loc>() ? 1 : 0 ) };

	/* horizontal bottom edges for top blocks -> MUSCL */
	if constexpr ( atTop<loc>() ){
	for ( Index x=0; x<block.nx(); ++x ){
		HLLC(
			RT::secondOrder, ER::horizontal,
			block.hBeg() + x,
			block.template lrBeg_block<ER::horizontal>() + x
		);
	}}
	/* inner horizontal edges */
	for ( Index y=y0; y<ny+1      ; ++y ){
	for ( Index x= 0; x<block.nx(); ++x ){
		HLLC(
			( recType<loc, ER::horizontal>() ), ER::horizontal,
			block.hBeg() + x + y*block.ih().mesh.nx,
			block.template lrBeg_block<ER::horizontal>() + x + y*block.nxRead()
		);
	}}
	/* horizontal top edges for bottom blocks -> MUSCL */
	if constexpr ( atBottom<loc>() ){
		Index
			nEdgesToAdd { block.ny() * block.ih().mesh.nx },
			nCellsToAdd { block.ny() * block.nxRead() };

		for ( Index x=0; x<block.nx(); ++x ){
			HLLC(
				RT::secondOrder, ER::horizontal,
				block.hBeg() + x + nEdgesToAdd,
				block.template lrBeg_block<ER::horizontal>() + x + nCellsToAdd
			);
		}
	}

	for ( Index y=0; y<block.ny(); ++y ){
		/* vertical left edges for right blocks -> MUSCL */
		if constexpr ( atRight<loc>() ){
			HLLC(
				RT::secondOrder, ER::vertical,
				block.vBeg() + y*(block.ih().mesh.nx + 1),
				block.template lrBeg_block<ER::vertical>() + y*block.nxRead()
			);
		}
		/* inner vertical edges */
		for ( Index x=x0; x<nx+1; ++x ){
			HLLC(
				( recType<loc, ER::vertical>() ), ER::vertical,
				block.vBeg() + x + y*(block.ih().mesh.nx + 1),
				block.template lrBeg_block<ER::vertical>() + x + y*block.nxRead()
			);
		}
		/* vertical right edges for left blocks -> MUSCL */
		if constexpr ( atLeft<loc>() ){
			HLLC(
				RT::secondOrder, ER::vertical,
				block.vBeg() + block.nx() + y*(block.ih().mesh.nx + 1),
				block.template lrBeg_block<ER::vertical>() + block.nx() + y*block.nxRead()
			);
		}
	}
	// std::cout << s.str() << std::endl;

	#undef HLLC
}

} // namespace hms::swe

#include "bi_computeFluxes_vec.tpp"
#include "bi_computeFluxes_cell.tpp"
#include "bi_fluxes.tpp"
#include "bi_fluxes_vec.tpp"
